Apparatus and method for estimating attitude using inertial measurement equipment

ABSTRACT

An apparatus for estimating attitude using an inertia measuring apparatus comprising: a Kalman filter  3  adapted to receive a quantity of state ω(t) from the inertia measuring apparatus  1 , input the quantity of state into a predetermined system equation and a predetermined observation equation and execute time update processing and observation update processing with regard to the system equation and the observation equation; a processing section  4  adapted to generate an estimated value of a time differentiation of a modified Rodrigues parameter α(t) based upon an output of the Kalman filter; an integral processing section  5  adapted to generate an estimated value of the modified Rodrigues parameter based upon an output of the processing section; a system propagation matrix generating section  6  adapted to update an observation sensitivity matrix based upon an initial value being given beforehand and an output of the integral processing section and supply the updated output to the Kalman filter; a transformed matrix generating section  7  adapted to generate a coordinate transformation matrix R based upon an output of the integral processing section; and an attitude estimation section  8  adapted to carry out an attitude estimation based upon an output of the transformed matrix generating section, is disclosed.

CROSS REFERENCE TO RELATED APPLICATION

This application is a 35 U.S.C. §371 of and claims priority to PCTInternational Application Number PCT/JP02/04103, which was filed 24 Apr.2002 (24.04.02), and was published in Japanese which was based onJapanese Patent Application No. 2001-128261 which was filed 24 Apr. 2001(25.04.01) and the teachings of which are incorporated herein byreference.

FIELD OF THE INVENTION

The present invention relates to an apparatus for estimating attitude,which, for example, estimates the attitude of an artificial satellitebased upon signals from an inertia measuring apparatus, and to a methodfor estimating the attitude and a program for the same.

BACKGROUND ART

In recent years, in the field of space development, attentions have beenpaid to a microsat (micro-miniaturized artificial satellite), which ismonofunctional but can be developed easily, in addition to big andmulti-functional artificial satellites. It is considered that easinessin the development of such microsat can reduce the cost per unit of anartificial satellite and also renders the space development morefamiliar. At present, a great number of microsats have been launched aspiggyback satellites.

The merits of the microsat are the features of monofunctional and lowcost. However, in order to provide the microsat with high performance,it is essential for the satellite to be equipped with three axialcontrol mechanism. Because of that, an attitude observation system withhigh performance is required. However, considering the increase in thedevelopmental cost for the satellite, it is not possible to use such aexpensive hardware.

Conventionally, various types of sensors have been used as an attitudesensor for the three axial satellites. For example, it has been knownthat exact attitude estimation including rolls and yaws can be done bydetecting an attitude angular velocity from a gyro, such as an inertiareference apparatus, and then integrating the attitude angular velocity.The inertia reference apparatus is an apparatus which is mounted with arate gyro or the like in the three axial manner and packaged together.When directly integrating signals of the gyro and then using them, thereis a merit for obtaining dynamical and good responsive detection. Thenagain, there is a demerit of accumulation of drifts during a long run.In order to compensate such accumulation of drifts, a manner ofrendering feedback of the errors using a star tracker from the absolutereference to the detection system has been employed. However, becausethis manner requires expensive sensors, it is not possible to apply itto a microsat that is monofunctional and is required to be low cost.

The present invention was made for aiming at achieving theabove-mentioned object, and it is an object of the present invention toprovide an attitude estimation apparatus being provided with highaccuracy and inexpensive and a method for estimating attitude, i.e., anattitude estimation apparatus using an inertia measuring apparatus andbeing capable of estimating the attitude of an artificial satellite bymeans of hardware, such as a piezoelectric vibration gyro and a fibergyro, of which cost per unit is relatively low, and a method and aprogram for estimating the attitude.

With the present invention, it becomes possible to estimate the attitudeonly with information from the gyro without requiring sensors such as astar tracker. The conventional artificial satellites may face such acircumstance where they cannot use the sensors. For example, the presentinvention can be a countermeasure under a practical operation when anartificial satellite cannot use a sensor such as a star tracker in themaneuver.

DISCLOSURE OF THE INVENTION

The attitude estimation apparatus using an inertia measuring apparatusaccording to the present invention comprises:

a Kalman filter adapted to receive a quantity of state ω(t) from theinertia measuring apparatus to input it into a predetermined systemequation and observation equation and execute time update processing andobservation update processing with regard to the equations,

a processing section adapted to generate an estimated value of a timedifferentiation of a modified Rodrigues parameter α(t) based upon anoutput of the Kalman filter,

an integral processing section adapted to generate an estimated value ofa modified Rodrigues parameter based upon an output of the processingsection,

a system propagation matrix generating section adapted to update anobservation sensitivity matrix based upon an initial value givenbeforehand and an output of the integral processing section, and supplythe updated output to the Kalman filter,

a transformed matrix generating section adapted to generate a coordinatetransformation matrix R based upon an output of the integral processingsection, and

an attitude estimation section adapted to execute attitude estimationbased upon an output of the transformed matrix generating section.

The modified Rodrigues parameter is preferably given according toequations (2.2.31) and (2.2.32) described later wherein e is adirectional cosine of a rotative axis and θ is a rotation angle aboutthe rotative axis.

When an equation (2.3.1) described later is provided, the observationequation is preferably given by adding observation noises into anequation (2.3.2).

When the right member of an equation (2.3.5) described later is set tobe f(α(t), {dot over (α)}(t),t), the system equation is preferably givenby adding system noises into an equation (2.3.6) described later.

The method for estimating attitude according to the present invention,wherein an inertia measuring apparatus is used, comprises:

first step for receiving signals from the inertia measuring apparatusand generating a system equation and an observation equation, wherein aquantity of state of a time differentiation of a modified Rodriguesparameter α(t) and a quantity of state ω(t) are quantities to beobserved,

second step for obtaining an estimated value of a time differentiationof the parameter α(t) by means of an extended Kalman filter,

third step for numerically integrating the estimated value obtained inthe second step and obtaining an estimated value of the parameter α(t),

fourth step for determining a transformed matrix R(t) based upon theestimated value obtained in the third step, and

a step for obtaining an estimated value of the attitude based upon thetransformed matrix R(t).

The program according to the present invention is to render a computerto execute:

first step for receiving signals from the inertia measuring apparatusand generating a system equation and an observation equation, in thosewhich a quantity of state of a time differentiation of a modifiedRodrigues parameter α(t) and a quantity of state ω(t) are quantities tobe observed,

second step for obtaining an estimated value of a time differentiationof the parameter α(t) by means of an extended Kalman filter,

third step for numerically integrating the estimated value obtained inthe second step to obtain an estimated value of the parameter α(t),

fourth step for determining the transformed matrix R(t) based upon theestimated value obtained in the third step, and

a step for obtaining an estimated value of the attitude based upon thetransformed matrix R(t).

The program according to the present invention is recorded in, forexample, a recording medium.

The recoding medium includes, for example, a flexible disk, a hard disk,a magnetic tape, a magnet-optical disk, a CD (including CD-ROM,Video-CD), a DVD (including DVD-Video, DVD-ROM, DVD-RAM), aROM-cartridge, a RAM memory cartridge equipped with battery backup, aflash memory cartridge, a nonvolatile RAM cartridge, and the like.

Furthermore, the recording medium also includes communication media,such as a cable communication medium, e.g., telephone lines, and awireless communication medium, e.g., microwave lines. Internet is alsoincluded in the communication media defined herein.

The medium is defined as an object in which information (mainly digitaldata and programs) are recorded by means of certain physical means andcapable of rendering a processing device, such as a computer and adedicated processor, to execute a preset function. In short, anymaterial can be the medium as far as it can execute downloading of aprogram into a computer by any means and it can render the computer toimplement the program.

The present invention is constituted by applying the modified Rodriguesparameter to the extended Kalman filter.

In general, a conversion equation using quaternions is known as apropagation equation of an attitude. In case of the quaternions, thereare such advantages that no specific point with respect to the attitudeexists and no transcendental function is used, which are unlike theEuler angles. Thus, the quaternions are generally used for describingthe attitude of a satellite. However, when the quaternions are estimatedwith respect to information only from a gyro, a required observationequation is not established in the Kalman filter that is an estimationalgorithm. In order to carry out the estimation by means of the Kalmanfilter, it is required to establish required equations, i.e., a systemequation and an observation equation. Hence, in the present invention,it has become possible to establish those equations by using themodified Rodrigues parameter that is an application of the quaternions.

BRIEF EXPLANATION OF THE DRAWINGS

FIG. 1 is a block diagram showing the apparatus according to oneembodiment for carrying out the present invention.

FIG. 2 is a flow chart showing the processing according to oneembodiment for carrying out the present invention.

FIG. 3 shows an example of a result of simulation.

BEST MODE FOR CARRYING OUT THE INVENTION

Now, the apparatus and the method according to embodiments for carryingout the present invention are explained with reference to the appendeddrawings. The explanation in the following is an example where theapparatus and the method according to embodiments for the presentinvention are used for estimating an attitude. The apparatus and themethod according to embodiments for the present invention are intendedto carry out estimation of the attitude information of an artificialsatellite for long run with high accuracy based upon only an output froma gyro serving as a sensor. In the present invention, modified Rodriguesparameters are used as the attitude determination equation of anartificial satellite.

FIG. 1 is a block diagram showing the apparatus according to oneembodiment for the present invention. The apparatus receives a quantityof state (angular velocity) ω(t) from a rate gyro serving as an inertiameasuring apparatus, and Kalman filter 3 generates a system equation andan observation equation and executes time update processing andobservation update processing according to the following equations.Besides, the details of the system equation and the observation equationwill be explained later.

Time Update{umlaut over ( α=f({dot over ( α(t),t){dot over ( α(t _(k−1))={dot over({circumflex over (α)}(t _(k−1)){dot over (P)} (t)=F(t) P (t)+ P (t)F(t)^(T) +B(t)Q(t)B(t)Observation Update{dot over ({circumflex over (α)}(t _(k))={dot over ( α(t _(k))+K[ω(t_(k))−h({dot over ( α(t _(k)))]K(t _(k))= P (t _(k))H(t _(k))^(T) [H(t _(k)) P (t _(k))H(t _(k))^(T)+R(t _(k))]⁻¹P(t _(k))₊ =[I−K(t _(k))H(t _(k))] P (t _(k))

In the above equations, a mark + shown with a subscript letterrepresents a variable after the observation update processing carriedout by the Kalman filter. Similarly, k shown with a subscript letterrepresents that implementation of an expression which is discretelyprocessed at an interval Δt is the k-th in the order, and k+1 representsthat the implementation is the k+1-th in the order. A negative numeralof −1 shown with a superscript letter represents an inverse matrix. Tshown with a superscript letter represents a transposed matrix.

Now, brief explanation will be made about the time update and theobservation update. The time update is defined as a processing tointegrate a linear model expressed in the following equation.{dot over (x)}(t)=F(t)x(t)+B(t)u(t)In the above-described case, the time update is a processing tointegrate the following equation.{umlaut over (α)}=F(t){dot over (α)}+B(t)u(t)This is a processing to obtain a predicted value by conductingcalculation. In addition, a covariance matrix of errors issimultaneously updated. The equation with regard to the covariancematrix of errors is expressed as the following equation.{dot over (P)} (t)=F(t) P (t)+ P (t)F(t)^(T) +B(t)Q(t)B(t)^(T)When conducting the update processing, this equation is subjected tointegral processing. As a result, P(t) is worked out.

The observation update is defined as a processing to work out a newestimated value based on an observed value and a predicted value byconducting calculation, those which are newly obtained. In theabove-described case, the equation is expressed as follows.{dot over ({circumflex over (α)}(t _(k))₊={dot over ( α(t _(k))+K(t_(k))[ω(t _(k))−H(t _(k)){dot over ( α(t _(k))]The component attached with a hat represents an estimated value, and thecomponent with a bar represents a predicted value. In this processing,information with regard to errors is obtained by comparing an observedvalue and a predicted value, the information is then multiplied by acertain gain to add the multiplied information to the predicted value,and the predicted value having been added is used as the estimatedvalue. Particularly, the part ofω(t_(k))−H(t_(k)){dot over ( α(t_(k))is called as innovation process or O-C, which is a process to comparethe observed value with the observed value predicted by conductingcalculation. With the process, a comparison is made to see if thepredicted value by conducting calculation and the observed value differwith each other or not. Information as regard to errors is obtained fromthe comparison, however, the information includes both usefulinformation and meaningless information. So-called Kalman gain expressedwithK(t_(k))is for determining how to divide such information and add the dividedinformation to the predicted value.

The Kalman gain is worked out according to the following equation.K(t _(k))= P (t _(k))H(t _(k))^(T) [H(t _(k)) P (t _(k))H(t _(k))^(T)+R(t _(k))]⁻¹(For conducting calculation according to the above equation, it wasrequired to determine P(t) by means of the time update.) Besides,R(t_(k)) is one called as a covariance matrix of the observation noises.

Further, the covariance matrix of errors is also updated with the Kalmangain, the update is carried out according to the following equation.P(t _(k))₊ =[I−K(t _(k))H(t _(k))] P (t _(k))

The “time update” and “observation update” those which execute theprocessings as described above are used to estimate the quantity ofstate, such as “position” and “attitude”, by the time update during thetime after executing the observation update with a new observed valueuntil obtaining a new observed value.

A processing section 4 generates {dot over ({circumflex over (α)} (anestimated value of a time differentiation {dot over (α)}(t) of amodified Rodrigues parameter α(t)) based upon an output of the Kalmanfilter 3.

An integral processing section 5 generates estimated values {circumflexover (α)}, {circumflex over (α)}₀ of the modified Rodrigues parameterbased upon an output of the processing section 4.

A system propagation matrix generating section 6 updates an observationsensitivity matrix based upon an output of an initial value settingsection 2 and an output of the integral processing section 5 andsupplies the updated output to the Kalman filter 3. The systempropagation matrix is F(t) in a linear system expressed as an equationas follows.{dot over (x)}(t)=F(t)x(t)+B(t)u(t) (In this case:{umlaut over(α)}=F(t){dot over (α)}+B(t)u(t))In this case, the system propagation matrix is obtained by subjectingthe quantity of state {dot over (α)} of an equation (2.3.5) describedlater to a partial differentiation (equations (2.3.8) to (2.3.16)described later correspond to the components). The system propagationmatrix is a matrix representing a relation between the quantity of stateand a differentiated value of the quantity of state. The observationsensitivity matrix is H(t) expressed in the following equation.ω(t)=H(t){dot over (α)}(t)+v(t)

The observation sensitivity matrix is a matrix representing a relationbetween an observed value and the quantity of state (The matrix isobtained from an equation (2.3.1) described later. However, since thisvalue cannot be solved holomorphically, it is subjected to a numericaloperation processing).

The transformed matrix generating section 7 generates a coordinatetransformation matrix R based on an output of the integral processingsection 5.

The attitude estimation section 8 executes attitude estimation basedupon an output of the transformed matrix generating section 7.

Now, the operation will be explained. The outline of the processings inthe present invention will be explained first, the theoretics and thedetails of algorithm of the present invention will be then explained,and the result of simulation will be finally explained.

In an embodiment for carrying out the present invention, a knownextended Kalman filter is used. The Kalman filter is a manner where aninitial value of a desired quantity of state is previously given andestimated values of the quantity of state are sequentially calculatedbased upon the initial value and the time record of the quantity ofobservation.

FIG. 2 shows an outlined flow chart of the operation.

-   S1: Preparing a system equation and an observation equation by    receiving signals from a rate gyro 1, in those which a quantity of    state of a time differentiation {dot over (α)}(t) of a modified    Rodrigues parameter α(t) and ω(t) are the quantities to be observed.-   S2: Working out an estimated value of {dot over (α)}(t) with the    extended Kalman filter.-   S3: Working out an estimated value of α(t) by numerically    integrating the estimated value obtained in the step S2.-   S4: Determining the transformed matrix R(t) based upon the estimated    value obtained in the step S3.-   S5: Working out an estimated value of the attitude at the time t    based upon the transformed matrix R(t).

As described above, the system equation and the observation equation canbe generated by using the modified Rodrigues parameter, and it becamepossible to apply the extended Kalman filter for the attitude estimationof an artificial satellite. As a result, the influence caused by driftswas reduced, and it became possible to carry out the attitude estimationwith high accuracy by means of only the rate gyro. Inversely, it was notpossible to generate the system equation and the observation equation bythe conventional manners, and the extended Kalman filter could not beapplied to the attitude estimation of an artificial satellite.

Now, explanation will be made on the theoretics.

[Theoretics]

1. Rotation of Vectors About an Arbitrary Axis

Now, assuming that a certain vector is rotated through angle θ about anarbitrary axis. Here, a unit vector in the axial direction is designatedby e, and the vector after the rotation is designated by x′. In thatcase, the following equation becomes valid.x′=(x·e)e+{x−(x·e)e} cos θ+(e×x)sin θ  (2.1.1a)=(x·e)e−x+x+{x−(x·e)e} cos θ+sin θ(e×x)=x+(1−cos θ){(x·e)e−x}+sin θ(e×x)=x+(1−cos θ){(x·e)−(e·e)x}+sin θ(e×x)=x+(1−cos θ)e×e×x+sin θ(e×x)={I+sin θê+(1−cos θ)ê ² }x  (2.1.1b)wherein I represents a unit matrix. Further, under the followingconditions,

$\begin{matrix}{\hat{e} = \left\lbrack {\begin{matrix}0 \\e_{z} \\{- e_{y}}\end{matrix}\begin{matrix}{- e_{z}} \\0 \\e_{x}\end{matrix}\begin{matrix}e_{y} \\{- e_{x}} \\0\end{matrix}} \right\rbrack} & \left( {2.1{.2}a} \right)\end{matrix}$êx=e×x  (2.1.2b)

the following equation becomes valid based on the formulae of doubleangles and half angles.

$\begin{matrix}{x^{\prime} = {\left\{ {I + {2\cos\;\frac{\theta}{2}\sin\;\frac{\theta}{2}\hat{e}} + {2\sin\;\frac{\theta}{2}{\hat{e} \cdot \sin}\;\frac{\theta}{2}\hat{e}}} \right\} x}} & \left( {2.1{.3}} \right)\end{matrix}$

Here, if the quaternion q that represents the rotation of an arbitraryaxial rotation θ is expressed as follows,

$\begin{matrix}{q = \left( {{\cos\;\frac{\theta}{2}},{\sin\;\frac{\theta}{2}e}} \right)} & \left( {2.1{.4}a} \right)\end{matrix}$=(α,V)  (2.1.4b)

x′ is expressed as follows.x′={I+2α{circumflex over (V)}+2{circumflex over (V)} ² }x=Rx  (2.1.5)However, in the above equation, the following equations are valid.

$\begin{matrix}{{\hat{V} = \left\lbrack {\begin{matrix}0 \\v_{z} \\{- v_{y}}\end{matrix}\begin{matrix}{- v_{z}} \\0 \\v_{x}\end{matrix}\begin{matrix}v_{y} \\{- v_{x}} \\0\end{matrix}} \right\rbrack}{R = \left\{ {I + {2a\hat{V}} + {2{\hat{V}}^{2}}} \right\}}} & \left( {2.1{.6}} \right)\end{matrix}$

Hence, if the direction of the axis at the time to and the transformedmatrix R at the time t are known, the direction of the axis at the timet, namely the attitude, can be known.

2. System Equation

If basis vectors in tridimensional space are c₁, C₂, C₃, and when thesethree vectors rotate about the axis that is expressed with a certainangular velocity vector shown in the following,

$\begin{matrix}{\omega = \begin{bmatrix}\omega_{1} \\\omega_{2} \\\omega_{3}\end{bmatrix}} & \left( {2.2{.1}} \right)\end{matrix}$ω=ω₁ c ₁+ω₂ c ₂+ω₃ c ₃  (2.2.2)

the following relation comes to valid with respect to the timeintegrations of c₁, c₂, and c₃.

$\begin{matrix}\left\{ \begin{matrix}{{\overset{.}{c}}_{1} = {\omega \times c_{1}}} \\{{\overset{.}{c}}_{2} = {\omega \times c_{2}}} \\{{\overset{.}{c}}_{3} = {\omega \times c_{3}}}\end{matrix} \right. & \left( {2.2{.3}} \right)\end{matrix}$

When the equation (2.2.2) is collected into one equation by using matrixexpression and using {circumflex over (ω)} indicated as follows,

$\begin{matrix}{\hat{\omega} = \left\lbrack {\begin{matrix}0 \\\omega_{3} \\{- \omega_{2}}\end{matrix}\begin{matrix}{- \omega_{3}} \\0 \\\omega_{1}\end{matrix}\begin{matrix}\omega_{2} \\{- \omega_{1}} \\0\end{matrix}} \right\rbrack} & \left( {2.2{.4}} \right)\end{matrix}$the equation comes to an equation as follows.ċ=c{circumflex over (ω)}  (2.2.5)wherein,c=[c ₁ |c ₂ |c ₃]  (2.2.6)c{circumflex over (ω)}=[ω×c ₁ ω×c ₂ |ω×c ₃]  (2.2.7)3. Use of Modified Rodrigues Parameter and Its First-Order DifferentialEquation

Modified Rodrigues parameter: α and α₀ as shown below are defined asfollows.

$\begin{matrix}{{a \equiv {4\tan\;\frac{\theta}{4}e}} = \begin{bmatrix}\alpha_{x} \\\alpha_{y} \\\alpha_{z}\end{bmatrix}} & \left( {2.2{.31}} \right) \\{\alpha_{0} \equiv {2 - \frac{a^{2}}{8}}} & \left( {2.2{.32}} \right)\end{matrix}$

That is, α is a tridimensional vector. e represents a directional cosineof the rotative axis, and θ represents a rotation angle about therotative axis. If the modified Rodrigues parameter defined according tothe equation (2.2.31) is used, the element of the quaternion in theequation (2.1.4) can be expressed with α. If a and V expressed with αare used, the transformed matrix of the equation (2.1.5) can beexpressed with α.

Further, the time history of R can be worked out by preparing thefirst-order differential equation of α. The first-order differentialequation of α is merely three simultaneous equations.

From the equation (2.2.32), a is worked out as follows.

$\begin{matrix}\begin{matrix}{a^{2} = {16 - {8\alpha_{0}}}} \\{a = {{\cos\;\frac{\theta}{2}} = {{2\cos^{2}\frac{\theta}{4}} - 1}}} \\{= {\frac{2}{1 + {\tan^{2}\frac{\theta}{4}}} - 1}} \\{= {{\frac{2}{1 + \frac{a^{2}}{16}} - 1} = {{\frac{32}{16 + a^{2}} - 1} = {\frac{32}{16 + 16 - {8\alpha_{0}}} - 1}}}} \\{= {\frac{4}{4 - \alpha_{0}} - 1}}\end{matrix} & \left( {2.2{.33}} \right)\end{matrix}$Hence, the following equation comes valid.

$\begin{matrix}{a = \frac{\alpha_{0}}{4 - \alpha_{0}}} & \left( {2.2{.34}} \right)\end{matrix}$

Next, attention is paid to V. V is worked out according to the followingequation.

$\begin{matrix}{V = {{\sin\;\frac{\theta}{2}e} = {{\frac{4\tan\;\frac{\theta}{4}}{4 - \alpha_{0}}e} = \frac{a}{4 - \alpha_{0}}}}} & \left( {2.2{.35}} \right)\end{matrix}$

The transformed matrix R is worked out according to the followingequations based on the equations (2.2.34) and (2.2.35).

$\begin{matrix}\begin{matrix}{R = {I + {2a\hat{V}} + {2{\hat{V}}^{2}}}} \\{= {I + {\frac{2\alpha_{0}}{4 - \alpha_{0}} \cdot \frac{\hat{a}}{4 - \alpha_{0}}} + \frac{2{\hat{a}}^{2}}{\left( {4 - \alpha_{0}} \right)^{2}}}} \\{= {I + {\frac{2\alpha_{0}}{\left( {4 - \alpha_{0}} \right)^{2}}\hat{a}} + {\frac{2}{\left( {4 - \alpha_{0}} \right)^{2}}{\hat{a}}^{2}}}}\end{matrix} & \left( {2.2{.36}} \right)\end{matrix}$

If working out the first-order differential equation of a using theequations (2.2.34) and (2.2.35), the following equation comes valid.

$\begin{matrix}\begin{matrix}{\overset{.}{a} = {\left( {{\frac{\alpha_{0}}{2}I} + {\frac{1}{2}\hat{a}} + {\frac{1}{8}{aa}^{T}}} \right)\omega}} \\{= {\frac{1}{8}\left\{ {\left\lbrack {\begin{matrix}{4\alpha_{0}} \\0 \\0\end{matrix}\begin{matrix}0 \\{4\alpha_{0}} \\0\end{matrix}\begin{matrix}0 \\0 \\{4\alpha_{0}}\end{matrix}} \right\rbrack + \left\lbrack {\begin{matrix}0 \\{4\alpha_{z}} \\{{- 4}\alpha_{2}}\end{matrix}\begin{matrix}{- \alpha_{z}} \\0 \\{4\alpha_{x}}\end{matrix}\begin{matrix}\alpha_{y} \\{- \alpha_{x}} \\0\end{matrix}} \right\rbrack +} \right.}} \\{\begin{bmatrix}\alpha_{x} \\\alpha_{y} \\\alpha_{z}\end{bmatrix}\left\lbrack {\begin{matrix}\alpha_{x}\end{matrix}\begin{matrix}\alpha_{y}\end{matrix}\begin{matrix}{\left. \left. \alpha_{z} \right\rbrack \right\}\omega}\end{matrix}} \right.} \\{= {\frac{1}{8}\left\{ {\begin{bmatrix}{4\alpha_{0}} & {{- 4}\alpha_{z}} & {4\alpha_{y}} \\{4\alpha_{z}} & {4\alpha_{0}} & {{- 4}\alpha_{x}} \\{{- 4}\alpha_{y}} & {4\alpha_{x}} & {4\alpha_{0}}\end{bmatrix} + \begin{bmatrix}\alpha_{x}^{2} & {\alpha_{x}\alpha_{y}} & {\alpha_{x}\alpha_{z}} \\{\alpha_{x}\alpha_{y}} & \alpha_{y}^{2} & {\alpha_{y}\alpha_{z}} \\{\alpha_{x}\alpha_{z}} & {\alpha_{y}\alpha_{z}} & \alpha_{z}^{2}\end{bmatrix}} \right\}\omega}} \\{{\therefore\overset{.}{a}} = {{\frac{1}{8}\begin{bmatrix}{{4\alpha_{0}} + \alpha_{x}^{2}} & {{{- 4}\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\alpha_{y}} + {\alpha_{x}\alpha_{z}}} \\{{4\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\alpha_{0}} + \alpha_{y}^{2}} & {{{- 4}\alpha_{x}} + {\alpha_{y}\alpha_{z}}} \\{{{- 4}\alpha_{y}} + {\alpha_{x}\alpha_{z}}} & {{4\alpha_{x}} + {\alpha_{y}\alpha_{z}}} & {{4\alpha_{0}} + \alpha_{z}^{2}}\end{bmatrix}}\omega}}\end{matrix} & \left( {2.2{.41}} \right)\end{matrix}$[Details of Algorithm]4. Application to Extended Kalman Filter

In an embodiment for the present invention, the attitude estimation iscarried out based on the angular velocity vector ω(t) measured from therate gyro which is mounted to three axes being perpendicular to oneanother in a satellite.

More particularly, the quantity of state of the time differentiation{dot over (α)}(t) of a modified Rodrigues parameter α(t), and the systemequation and the observation equation in those which ω(t) is thequantity of observation are prepared from the equation (2.2.41), and theestimated value of {dot over (α)}(t) is worked out with the extendedKalman filter. Then, the obtained estimated value is numericallyintegrated to obtain the estimated value of α(t), which is thensubstituted into the equation (2.2.36) to determine the transformedmatrix R(t), and the estimated value of the attitude at the time t isworked out according to the equation (2.1.5).

In this case, x in the right member of the equation (2.1.5) is theinitial value of the directional vector that represents an arbitraryfixed axis of an artificial satellite, and x′ in the left member becomesthe direction of the fixed axis after t seconds. The norms of x′ and xare equivalent.

In the first place, the observation equation is worked out. If set asfollows based on the equation (2.2.41),

$\begin{matrix}{H^{- 1} = \begin{bmatrix}{{4\alpha_{0}} + \alpha_{x}^{2}} & {{{- 4}\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\alpha_{y}} + {\alpha_{x}\alpha_{z}}} \\{{4\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\alpha_{0}} + \alpha_{y}^{2}} & {{{- 4}\alpha_{x}} + {\alpha_{y}\alpha_{z}}} \\{{{- 4}\alpha_{y}} + {\alpha_{x}\alpha_{z}}} & {{4\alpha_{x}} + {\alpha_{y}\alpha_{z}}} & {{4\alpha_{0}} + \alpha_{z}^{2}}\end{bmatrix}} & \left( {2.3{.1}} \right)\end{matrix}$the following equation becomes valid.ω=H{dot over (α)}  (2.3.2)

By adding the observation noises to the equation (2.3.2), theobservation equation in the extended Kalman filter is obtained.

Now, the system equation will be worked out. The following equation isobtained from the equation (2.2.32).

$\begin{matrix}{\alpha_{0} = {{{- \frac{1}{4}}\left( {a \cdot \overset{.}{a}} \right)} = {{- \frac{1}{4}}\left( {{\alpha_{x}{\overset{.}{\alpha}}_{x}} + {\alpha_{y}{\overset{.}{\alpha}}_{y}} + {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)}}} & \left( {2.3{.3}} \right)\end{matrix}$

When developing the equation (2.2.41), the following equation comesvalid.

$\begin{matrix}{\overset{.}{a} = {\frac{1}{8}\begin{bmatrix}{{\left( {{4\alpha_{0}} + \alpha_{x}^{2}} \right)\omega_{x}} + {\left( {{{- 4}\alpha_{z}} + {\alpha_{x}\alpha_{y}}} \right)\omega_{y}} + {\left( {{4\alpha_{y}} + {\alpha_{x}\alpha_{z}}} \right)\omega_{z}}} \\{{\left( {{4\alpha_{z}} + {\alpha_{x}\alpha_{y}}} \right)\omega_{x}} + {\left( {{4\alpha_{0}} + \alpha_{y}^{2}} \right)\omega_{y}} + {\left( {{{- 4}\alpha_{x}} + {\alpha_{y}\alpha_{z}}} \right)\omega_{z}}} \\{{\left( {{{- 4}\alpha_{y}} + {\alpha_{x}\alpha_{z}}} \right)\omega_{x}} + {\left( {{4\alpha_{x}} + {\alpha_{y}\alpha_{z}}} \right)\omega_{y}} + {\left( {{4\alpha_{0}} + \alpha_{z}^{2}} \right)\omega_{z}}}\end{bmatrix}}} & \left( {2.3{.4}} \right)\end{matrix}$If considering that the time differentiation dω/dt of ω is meaninglessand deleting it, the first-order differentiation of the equation (2.3.4)comes to as follows.

$\begin{matrix}{\overset{¨}{a} = {\frac{1}{8}\begin{bmatrix}{{\left( {{\alpha_{x}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{{- 4}{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{y}} + {\left( {{4{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}} \\{{\left( {{4{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{z}} + {\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} + {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} + {\left( {{{- 4}{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}} \\{{\left( {{{- 4}{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{4{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} + {\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} + {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}}\end{bmatrix}}} & \left( {2.3{.5}} \right)\end{matrix}$If the right member of the above equation is f(α(t),{dot over(α)}(t),t), the following equation comes valid.

$\begin{matrix}{{\overset{¨}{a}(t)} = {{f\left( {{a(t)},{\overset{.}{a}(t)},t} \right)} = \begin{bmatrix}f_{x} \\f_{y} \\f_{z}\end{bmatrix}}} & \left( {2.3{.6}} \right)\end{matrix}$By adding the system noises to the equation (2.3.6), the equationbecomes a system equation in the extended Kalman filter.

Since the equation (2.3.6) is a nonlinear differential equation, it isrequired to work out a matrix F which is obtained by partiallydifferentiating the equation (2.3.5) by the quantity-of-state vector inorder to work out the system propagation matrix.

$\begin{matrix}{F = {\left\lbrack \frac{\partial f}{\partial\overset{.}{a}} \right\rbrack = \begin{bmatrix}{{\partial f_{x}}/{\partial{\overset{.}{\alpha}}_{x}}} & {{\partial f_{x}}/{\partial{\overset{.}{\alpha}}_{y}}} & {{\partial f_{x}}/{\partial{\overset{.}{\alpha}}_{z}}} \\{{\partial f_{y}}/{\partial{\overset{.}{\alpha}}_{x}}} & {{\partial f_{y}}/{\partial{\overset{.}{\alpha}}_{y}}} & {{\partial f_{y}}/{\partial{\overset{.}{\alpha}}_{z}}} \\{{\partial f_{z}}/{\partial{\overset{.}{\alpha}}_{x}}} & {{\partial f_{z}}/{\partial{\overset{.}{\alpha}}_{y}}} & {{\partial f_{z}}/{\partial{\overset{.}{\alpha}}_{z}}}\end{bmatrix}}} & \left( {2.3{.7}} \right)\end{matrix}$

From the equation (2.3.5), each of the elements of F is obtained asfollows.

$\begin{matrix}{\frac{\partial f_{x}}{\partial{\overset{.}{\alpha}}_{x}} = {\left( {{\alpha_{x}\omega_{x}} + {\alpha_{y}\omega_{y}} + {\alpha_{z}\omega_{z}}} \right)/8}} & \left( {2.3{.8}} \right) \\{\frac{\partial f_{x}}{\partial{\overset{.}{\alpha}}_{y}} = {\left( {{{- \alpha_{y}}\omega_{x}} + {\alpha_{x}\omega_{y}} + {4\omega_{z}}} \right)/8}} & \left( {2.3{.9}} \right) \\{\frac{\partial f_{x}}{\partial{\overset{.}{\alpha}}_{z}} = {\left( {{{- \alpha_{z}}\omega_{x}} - {4\omega_{y}} + {\alpha_{x}\omega_{z}}} \right)/8}} & \left( {2.3{.10}} \right) \\{\frac{\partial f_{y}}{\partial{\overset{.}{\alpha}}_{x}} = {\left( {{\alpha_{y}\omega_{x}} - {\alpha_{x}\omega_{y}} - {4\omega_{z}}} \right)/8}} & \left( {2.3{.11}} \right) \\{\frac{\partial f_{y}}{\partial{\overset{.}{\alpha}}_{y}} = {\left( {{\alpha_{x}\omega_{x}} + {\alpha_{y}\omega_{y}} + {\alpha_{z}\omega_{z}}} \right)/8}} & \left( {2.3{.12}} \right) \\{\frac{\partial f_{y}}{\partial{\overset{.}{\alpha}}_{z}} = {\left( {{4\omega_{x}} - {\alpha_{z}\omega_{y}} + {\alpha_{y}\omega_{z}}} \right)/8}} & \left( {2.3{.13}} \right) \\{\frac{\partial f_{z}}{\partial{\overset{.}{\alpha}}_{x}} = {\left( {{\alpha_{z}\omega_{x}} + {4\omega_{y}} - {\alpha_{x}\omega_{z}}} \right)/8}} & \left( {2.3{.14}} \right) \\{\frac{\partial f_{z}}{\partial{\overset{.}{\alpha}}_{y}} = {\left( {{{- 4}\omega_{x}} + {\alpha_{z}\omega_{y}} - {\alpha_{y}\omega_{z}}} \right)/8}} & \left( {2.3{.15}} \right) \\{\frac{\partial f_{z}}{\partial{\overset{.}{\alpha}}_{z}} = {\left( {{\alpha_{x}\omega_{x}} + {\alpha_{y}\omega_{y}} + {\alpha_{z}\omega_{z}}} \right)/8}} & \left( {2.3{.16}} \right)\end{matrix}$[Simulation Result]

The result of simulation is shown in FIG. 3. FIG. 3 shows estimatederrors in the direction of the apparatus body axis when the standarddeviation of observation noises is equal to 3×10⁻⁵ [rad/sec]. Accordingto FIG. 3, the errors increase in proportion to time in the initialstage and reach 2.2 [deg] at 1.2×10⁵ [sec]. However, the errors are thenstabilized. According to the embodiment for the present invention, theKalman gain has converged soon. Hence, reconstruction of information canbe implemented in a good condition.

Moreover, as evident from FIG. 3, there is no accumulation of the driftsand the errors can be settled within a certain range in the embodimentfor the present invention. Contrary thereto, according to theconventional techniques, the drifts accumulate and the errors have neverbe settled within a certain range. In this concern, the technique of theembodiment for the present invention is excellently improved.

The present invention can be modified in various ways within a scopedescribed in the claims without being limited to the range of theembodiment described above. It is no need to say that such modificationsof the present invention also fall within the scope of the presentinvention.

The present invention can be applied not only to an artificial satellitebut also to the general fields that execute the estimation of attitudeby means of a gyro, that is any apparatuses those which use an inertiameasuring apparatus. The present invention can guarantee improvement inthe accuracy, and its application field will be over a wide range. Forexample, the present invention can be utilized for a function to preventa video camera shake and for attitude estimation of a car navigationsystem and the like.

Besides, in the specification for the present invention, it should benoted that the means described therein does not always represent aphysical means, and it also includes the functions of the respectivemeans those which are realized by software. Furthermore, the means mayfurther include such a case that the function of one means is realizedby two or more physical means, or the functions of two or more means arerealized by one physical means.

As described above, according to the present invention, the attitudeestimation with high accuracy can be realized based upon an output fromonly the inertia measuring apparatus. In particular, with the presentinvention, the drifts which have been problematic in the conventionaltechniques can be eliminated.

1. An apparatus for estimating attitude using an inertia measuringapparatus comprising: a Kalman filter adapted to receive a quantity ofstate ω(t) from the inertia measuring apparatus, input the quantity ofstate into a predetermined system equation and a predeterminedobservation equation, and execute time update processing and observationupdate processing with regard to the system equation and the observationequation, a processing section adapted to generate an estimated value ofa time differentiation of a modified Rodrigues parameter α(t) based uponan output of the Kalman filter an integral processing section adapted togenerate an estimated value of the modified Rodrigues parameter basedupon an output of the processing section, a system propagation matrixgenerating section adapted to update an observation sensitivity matrixbased upon an initial value being given beforehand and an output of theintegral processing section and supply the updated output to the Kalmanfilter, a transformed matrix generating section adapted to generate acoordinate transformation matrix R based upon an output of the integralprocessing section, and an attitude estimation section adapted to carryout an attitude estimation based upon an output of the transformedmatrix generating section.
 2. The apparatus for estimating attitudeusing an inertia measuring apparatus according to claim 1 characterizedin that the modified Rodrigues parameter is given according to thefollowing equations; $\begin{matrix}{{a \equiv {4\;\tan\;\frac{\theta}{4}e}} = \begin{bmatrix}\begin{matrix}\alpha_{x} \\\alpha_{y}\end{matrix} \\\alpha_{z}\end{bmatrix}} \\{\alpha_{0} \equiv {2 - \frac{a^{2}}{8}}}\end{matrix}$ wherein e is a directional cosine of a rotative axis and θis a rotation angle about the rotative axis.
 3. The apparatus forestimating attitude using an inertia measuring apparatus according toclaim 1 characterized in that the observation equation is given byadding observation noises to an equation;ω=H{dot over (α)} when the following equation is valid$H^{- 1} = {\begin{bmatrix}{{4\;\alpha_{0}} + \alpha_{x}^{2}} & {{{- 4}\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\;\alpha_{y}} + {\alpha_{x}\alpha_{z}}} \\{{4\;\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\;\alpha_{0}} + \alpha_{y}^{2}} & {{{- 4}\alpha_{x}} + {\alpha_{y}\alpha_{z}}} \\{{{- 4}\alpha_{y}} + {\alpha_{x}\alpha_{z}}} & {{4\;\alpha_{x}} + {\alpha_{y}\alpha_{z}}} & {{4\;\alpha_{0}} + \alpha_{z}^{2}}\end{bmatrix}.}$
 4. The apparatus for estimating attitude using aninertia measuring apparatus according to claim 1 characterized in thatthe system equation is given by adding system noises to the followingequation;${\overset{¨}{a}(t)} = {{f\left( {{a(t)},{\overset{.}{a}(t)},t} \right)} = \begin{bmatrix}f_{x} \\f_{y} \\f_{z}\end{bmatrix}}$ when the right member of the following equation;$\overset{¨}{a} = {\frac{1}{8}\begin{bmatrix}\begin{matrix}{{\left( {{\alpha_{x}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{{- 4}\;{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{y}} +} \\{\left( {{4{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}\end{matrix} \\\begin{matrix}{{\left( {{4\;{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{z}} + {\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} + {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} +} \\{\left( {{{- 4}\;{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}\end{matrix} \\\begin{matrix}{{\left( {{{- 4}\;{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{4\;{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} +} \\{\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}\end{matrix}\end{bmatrix}}$ is set to be f(a(t),{dot over (α)}(t),t).
 5. A methodfor estimating attitude using an inertia measuring apparatus comprising:first step adapted to receive signals from the inertia measuringapparatus and generate a system equation and an observation equation, inwhich a quantity of state of the time differentiation of the modifiedRodrigues parameter α(t) and the quantity of state ω(t) are thequantities of observation, second step adapted to work out an estimatedvalue of the time differentiation of the parameter α(t) with an extendedKalman filter, third step adapted to numerically integrate the estimatedvalue worked out in the second step and work out an estimated value ofthe parameter α(t), fourth step adapted to determine the transformedmatrix R(t) based upon the estimated value worked out in the third step,and a step adapted to work out an estimated value of the attitude basedupon the transformed matrix R(t).
 6. The method for estimating attitudeusing an inertia measuring apparatus according to claim 5 characterizedin that the modified Rodrigues parameter is given according to thefollowing equations; $\begin{matrix}{{a \equiv {4\;\tan\;\frac{\theta}{4}e}} = \begin{bmatrix}\begin{matrix}\alpha_{x} \\\alpha_{y}\end{matrix} \\\alpha_{z}\end{bmatrix}} \\{\alpha_{0} \equiv {2 - \frac{a^{2}}{8}}}\end{matrix}$ wherein e is a directional cosine of the rotative axis andθ is an rotation angle about the rotative axis.
 7. The method forestimating attitude using an inertia measuring apparatus according toclaim 5 characterized in that the observation equation is given byadding observation noises to the following equation;ω=H{dot over (α)} when the following equation is valid$H^{- 1} = {\begin{bmatrix}{{4\;\alpha_{0}} + \alpha_{x}^{2}} & {{{- 4}\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\;\alpha_{y}} + {\alpha_{x}\alpha_{z}}} \\{{4\;\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\;\alpha_{0}} + \alpha_{y}^{2}} & {{{- 4}\alpha_{x}} + {\alpha_{y}\alpha_{z}}} \\{{{- 4}\alpha_{y}} + {\alpha_{x}\alpha_{z}}} & {{4\;\alpha_{x}} + {\alpha_{y}\alpha_{z}}} & {{4\;\alpha_{0}} + \alpha_{z}^{2}}\end{bmatrix}.}$
 8. The method for estimating attitude using an inertiameasuring apparatus according to claim 5 characterized in that thesystem equation is given by adding system noises to the followingequation;${\overset{¨}{a}(t)} = {{f\left( {{a(t)},{\overset{.}{a}(t)},t} \right)} = \begin{bmatrix}f_{x} \\f_{y} \\f_{z}\end{bmatrix}}$ when the right member of the following equation;$\overset{¨}{a} = {\frac{1}{8}\begin{bmatrix}\begin{matrix}{{\left( {{\alpha_{x}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{{- 4}\;{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{y}} +} \\{\left( {{4{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}\end{matrix} \\\begin{matrix}{{\left( {{4\;{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{z}} + {\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} + {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} +} \\{\left( {{{- 4}\;{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}\end{matrix} \\\begin{matrix}{{\left( {{{- 4}\;{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{4\;{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} +} \\{\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}\end{matrix}\end{bmatrix}}$ is set to be f(a(t),{dot over (α)}(t),t).
 9. A computerreadable medium for storing a program rendering a computer to implement:first step adapted to receive signal from the inertia measuringapparatus and generate a system equation and an observation equation, inthose which the quantity of state of the time differentiation of themodified Rodrigues parameter α(t) and the quantity of state ω(t) are thequantities of observation, second step adapted to work out an estimatedvalue of the time differentiation of the parameter α(t) with theextended Kalman filter, third step adapted to numerically integrate theestimated value worked out in the second step and work out an estimatedvalue of the parameter α(t), fourth step adapted to determine thetransformed matrix R(t) based upon the estimated value worked out in thethird step, and a step adapted to work out an estimated value of theattitude based upon the transformed matrix R(t).
 10. The computerreadable medium according to claim 9 characterized in that the modifiedRodrigues parameter is given according to the following equations;$\begin{matrix}{{a \equiv {4\;\tan\;\frac{\theta}{4}e}} = \begin{bmatrix}\begin{matrix}\alpha_{x} \\\alpha_{y}\end{matrix} \\\alpha_{z}\end{bmatrix}} \\{\alpha_{0} \equiv {2 - \frac{a^{2}}{8}}}\end{matrix}$ wherein e is a directional cosine of the rotative axis andθ is a rotation angle about the rotative axis.
 11. The computer readablemedium according to claim 9 characterized in that the observationequation is given by adding observation noises to the followingequation;ω=H{dot over (α)} when the following equation is valid$H^{- 1} = {\begin{bmatrix}{{4\;\alpha_{0}} + \alpha_{x}^{2}} & {{{- 4}\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\;\alpha_{y}} + {\alpha_{x}\alpha_{z}}} \\{{4\;\alpha_{z}} + {\alpha_{x}\alpha_{y}}} & {{4\;\alpha_{0}} + \alpha_{y}^{2}} & {{{- 4}\alpha_{x}} + {\alpha_{y}\alpha_{z}}} \\{{{- 4}\alpha_{y}} + {\alpha_{x}\alpha_{z}}} & {{4\;\alpha_{x}} + {\alpha_{y}\alpha_{z}}} & {{4\;\alpha_{0}} + \alpha_{z}^{2}}\end{bmatrix}.}$
 12. The computer readable medium according to claim 9characterized in that the system equation is given by adding systemnoises to the following equation;${\overset{¨}{a}(t)} = {{f\left( {{a(t)},{\overset{.}{a}(t)},t} \right)} = \begin{bmatrix}f_{x} \\f_{y} \\f_{z}\end{bmatrix}}$ when the right member of the following equation;$\overset{¨}{a} = {\frac{1}{8}\begin{bmatrix}{{\left( {{\alpha_{x}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{{- 4}{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{y}} + {\left( {{4{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}} \\{{\left( {{4{\overset{.}{\alpha}}_{z}} + {{\overset{.}{\alpha}}_{x}\alpha_{y}} + {\alpha_{x}{\overset{.}{\alpha}}_{y}}} \right)\omega_{z}} + {\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} + {\alpha_{y}{\overset{.}{\alpha}}_{y}} - {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} + {\left( {{{- 4}{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}} \\{{\left( {{{- 4}{\overset{.}{\alpha}}_{y}} + {{\overset{.}{\alpha}}_{x}\alpha_{z}} + {\alpha_{x}{\overset{.}{\alpha}}_{z}}} \right)\omega_{x}} + {\left( {{4{\overset{.}{\alpha}}_{x}} + {{\overset{.}{\alpha}}_{y}\alpha_{z}} + {\alpha_{y}{\overset{.}{\alpha}}_{z}}} \right)\omega_{y}} + {\left( {{{- \alpha_{x}}{\overset{.}{\alpha}}_{x}} - {\alpha_{y}{\overset{.}{\alpha}}_{y}} + {\alpha_{z}{\overset{.}{\alpha}}_{z}}} \right)\omega_{z}}}\end{bmatrix}}$ is set to be f(a(t),{dot over (α)}(t),t).